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Abstract 

We analyse the thermal properties of a harmonic chain with weak 
correlated disorder. With the use of a perturbative approach we derive 
analytical expressions for the time-evolution of the chain temperature 
and of the heat flow when both ends of the chain are coupled to 
heat baths. Our analytical and numerical results demonstrate that 
specific long-range correlations of the isotopic disorder can suppress 
or enhance the vibrational modes in pre-defined frequency windows. 
In this way one can arrange a frequency-selective heat flow through 
disordered chains. 

Pacs: 44.10. +i, 63.50.Gh, 63.20.Pw 

1 Introduction 

One-dimensional (ID) harmonic chains constitute convenient models for the 
effective study of thermal transport in solids. On the one hand, they represent 
one of the few cases for which analytical results can be obtained with a 
relative ease. On the other hand, they can be used to describe the features of 



solids in the low-temperature region, where nonlinear effects can be neglected. 
For these reasons the heat conduction in ID harmonic chains has long been 
a subject of active investigation (see [H [2] and references therein). 

This paper is devoted to the study of the thermal properties of ID har- 
monic chains with correlated isotopic disorder. Since the mid-1990s the in- 
terest in models of this kind has been constantly increasing [3], especially 
after it became clear that specific long-range correlations in random poten- 
tials can create a sort of delocalisation transition even in ID models [H [5]. 
In the past ten years many studies, both at the theoretical [HI [71 [8] and 
experimental [U [10] level, have thrown light on the phenomenon of enhance- 
ment or suppression of the Anderson localisation in low-dimensional systems 
with correlated random potentials. Within the field of ID models with cor- 
related disorder, some recent studies have begun to analyse the chains with 
long-range correlated random masses. These early works were focused on the 
anomalous diffusion of an initially localised energy input [11], or on the struc- 
ture of the eigenmodes and on the statistics of the vibrational spectrum [T2] . 

The present work adopts a different perspective, considering harmonic 
chains as open systems weakly coupled at both ends to Langevin heat baths 
as in the Casher-Lebowitz model The purpose of our study is twofold: 
to analyse the approach to the steady-state regime of the chain and to dis- 
cuss the effects of long-range correlations of the isotopic disorder on the 
heat conduction. Accordingly, we first analyse the dynamics of ID harmonic 
chains with random masses, without specifying the statistical properties of 
the disorder. Using standard tools of stochastic calculus, we obtain ana- 
lytical expressions for the temperature profile and heat flux, describing the 
evolution of the system toward the steady-state regime. In the stationary 
limit, our expressions reduce to the formulae first obtained by Matsuda and 
Ishii [I^ and Lebowitz and coworkers [13j. 

After, we focus on the thermal effects of long-range correletions. In a 
homogeneous chain all vibrational modes are extended, whereas in a random 
chain with uncorrelated disorder only low-frequency modes are significantly 
non-localised. However, specific long-range correlations of the disorder can 
generate the delocalisation or enhanced localisation of vibrational modes in 
predefined frequency windows. Thus, with an appropriate selection of the 
frequency intervals one can create a random chain with a finite fraction of 
extended modes in the middle of the frequency spectrum. This allows one to 
have a random chain with specific properties of thermal transport peculiar 
to models that are neither completely random nor perfectly ordered. 
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As an example, we show how to construct a random chain whose conduc- 
tivity tends to scale with the size of the system regardless of the boundary 
conditions. We also demonstrate that with a proper choice of long-range cor- 
relations one can enhance the localisation of the eigenmodes within a wide 
range of the frequency spectrum. In this way one can reduce the thermal 
conductivity of the chain with respect to that of a chain with white-noise 
disorder. 



2 The temperature profile and heat flux 

We consider a harmonic chain of masses mi, . . . ,171^ and we study its 
dynamics when the chain is weakly coupled to two Langevin baths at tem- 
peratures T+ and T_. The time evolution is described by the system of 
stochastic Ito equations 
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+ 5ii^2\kBT_ dW_ + 5iN^2\kBT+ dW^ 



with i = 1, . . . ,N . In eq. ([T]), Pi and represent the momentum and the 
displacement from the equilibrium position of the z-th mass of the chain. 
The symbols W+{t) and W-{t) stand for two independent Wiener processes, 
while A represents the coupling constant of the chain to the baths. Finally, 
D is the tridiagonal force-constant matrix, with elements Djj = 2k5ij — 
kSij-i — kSij^i for i = 2, . . . , N — 1 and j = 1, . . . ,N while the first and 
last rows take different forms according to the chosen boundary conditions. 
For fixed boundary conditions one has JDij = 2k6ij — k6ij-i and Dnj = 
2k6N,j — k6N,j+i, while for free boundary conditions the matrix elements are 
Dij = kSij — kdij^i and 'Dnj = kd^j — kd^j+i- 

To analyse the dynamics of the chain it is convenient to switch to normal 
coordinates and momenta, defined by the equations 
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Here e^*"* is the k-th component of the i-th eigenvector of the rescaled force- 
constant matrix 

D= (M-i/2^)Tj^^_i/2 
with eigenvalue uf, i.e., De*^*) = ufe^^^ 



for the mass matrix M,- 



rriiSij 



In eq. (E]), the symbol M stands 
As can be seen from eq. ([2]), the k-th 
component of the i-th. vector e^*^ defines the amplitude of the oscillations of 
the fc-th mass in the i-th normal mode. For this reason we will refer to the 
vectors e^*^ as displacement eigenvectors. 

After passing to the normal coordinates, the dynamical equations of the 
chain can be written in a compact form by introducing the 2A^-component 
vector X = {Qi, . . . ,Qn, Pi, . . . , Pn) and the 2N x 2N matrix 



i,j-N 
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with the coupling matrix C being defined by 
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if hi = 1, 
otherwise 
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We also need to introduce the 2N x 2 matrix 
/ 
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^j2\kBT^/mie^i^ ^j2\kBTj^lmjqe 



(1) 

N 



(TV) 



(6) 



and the two-component Wiener process W(t) = iy+(t)). Then one 

can write the dynamical equation in the form 



ciX = AX + B <m. 



(7) 



Eqs. (jl]) and ([7]) show that the interaction of the chain with the bath couples 
the normal modes, which are not independent as would be the case for an 
autonomous chain. From eq. ([6]) one should also notice that the coupling of 
each normal mode to the baths is proportional to to the first and the A^-th 
components of the corresponding displacement vector. 
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Given the deterministic initial condition X(0) = Xq, the formal solution 
of the stochastic equation ([7]) can be expressed in the form 



X(t) = U(t) 



Xo+ / U(-r)Bc/W(r) 



(8) 



where U(t) = exp (At) is the fundamental matrix of eq. ([7]). Eq. (|8]) defines 
a stochastic process which is Markovian, Gaussian, and that becomes sta- 
tionary over sufficiently long time-scales provided that the eigenvalues of the 
matrix A have negative real parts [15]. In other words, the stochastic pro- 
cess dH]) is a Ornstein-Uhlenbeck process in the broad sense. The Gaussian 
character of the process implies that its statistical properties are completely 
determined once the first two moments are known. One can easily show 
that the average value has the form m(t) = (X(t)) = U(t)Xo , while the 
covariance matrix can be expressed as 

K,,(t) = ((X,(t)-m,(t))(X,(t)-m,(t))) 

Throughout this paper we use the symbol (■ ■ ■) to denote the average over 
realisations of the Wiener processes. 

According to eq. ([8]) the dynamics of the chain is completely determined 
by the evolution operator U(t). This operator can be evaluated perturba- 
tively, provided that the chain is weakly coupled to the baths, i.e., 

\uJi — UJj\ , , 

A « ' ' ^' (10) 

for i, j = 1, . . . , A^. From a physical point of view, the weak-coupling con- 
dition f lTOj) ensures the absence of resonance effects. Obviously, for eq. f lTO|) 
to be satisfied, the eigenfrequencies {ui} must be non- degenerate. Below we 
assume that the condition ( ITOi) is fulfilled. For this reason, we do not apply 
the techniques used by Visscher to analyse the case of degenerated frequency 
spectra [T6j . 

To obtain an expansion of U(t), we note that the latter is the solution of 
the matrix equation 

U = AU (11) 

with the initial condition U(0) = 1. Assuming that condition ( !T0|) is satisfied, 
one can split the matrix A in two terms 

A = Ao - AAi. (12) 
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Here the first term, i.e., [Aoj^^- = 5i,j~-N — '-^'j^i-Nj — OijAK^i^Nj- 
the generator of the unperturbed dynamics, modified by the correction term 
[Ai]j^. = {1 — 6ij) Ci-N,j-N ■ Physically, the exact decomposition f|T2l) means 
that one considers a set of N damped and independent oscillators as the 
unperturbed system, and regards the coupling among the oscillators as a 
perturbation. The formal solution of eq. flTT]) can now be obtained in terms 
of a time-ordered exponential. 



SijXCi 



-TV 



represents 



U(t) = exp (Aot) Texp 



-A 



(13) 



with A[^\t) = exp (— Aot) Ai exp (Aot) being the perturbative term in the 
interaction representation. Under the weak-coupling hypothesis (ITUl) . one 
can expand the time-ordered exponential in eq. f[T5]) and obtain 



U(t) = e^»* -X [ e^»(*-")Aie^°"rfr + O (A^) 
Jo 



(14) 



The knowledge of the correlation matrix enables one to compute the 
temperature profile of the chain, defined via the relation 



Pi 



N N 

EE 
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For the sake of simplicity, in what follows we restrict our attention to the 
case in which the masses of the chain are at rest in their equilibrium position 
at the initial time, i.e., we assume that Xq = 0. Inserting the momentum- 
momentum correlator corresponding to this initial condition in eq. fllSp . one 
eventually obtains that the temperature of the i-th mass of the chain is 
described by the following formula. 
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Here we have introduced the notations 7^ = for the damping coef- 

ficients and Qk 



'jI for the damping-shifted eigenfrequencies of the 



chain. Note that in the hmit t — )■ cxd eq. f|T6l) reproduces the Matsuda-Ishii 
formula for the stationary temperature profile. 

We now turn our attention to the heat flux. Applying the standard 
rules of stochastic calculus, one obtains that the average value of the energy 
ij(^) = 1/2 (P| + ulQl) associated to the fc-th vibrational mode evolves in 
time according to the equation 



dt 



where the incoming modal heat flow is given by 
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while the outgoing modal heat flow is 
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mi niM 

Note again that eqs. (ITS!) and (IT^ reduce to the Matsuda-Ishii formulae in 
the asymptotic limit t — t- 00. 



3 Correlated disorder and structure of nor- 
mal modes 

The above results were obtained for a chain composed of arbitrary masses 
mi, . . . ,mj^. We now focus our attention on the case of weakly random 
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masses, defined by the condition that the fluctuations Srrii of the masses 
around their common mean value M = fn^ be small, i.e., \J~{5rn^ <^ M. 

Here and below the symbol (■ ■ ■) denotes the average over realisations of the 
disorder. The weak- disorder assumption allows us to use a perturbative ap- 
proach and to work within the second-order approximation in the parameter 

[dniif' / M . In this scheme, it is enough to specify the statistical properties 



of the random succession {drrii] by giving the flrst two moments dnii = 



and {6miY = a"^ and by assuming that the normalised binary correlator 



X{1) = -==^ (20) 

is a known function. We assume that disorder is spatially invariant in the 
mean; therefore, the binary correlator (!20|) depends only on the index dif- 
ference /. We also suppose that x(0 is a decreasing and even function of 
I. 

When the chain is not coupled to the external baths, its dynamics obeys 
the equation rrinqn = DnfcQ'fc- Taking the Fourier transform of both sides 
with respect to time, this equation can be cast in the form 

Qn+l + qn-l + -j^^^n = ^2 - ^^^^ 9n- (21) 

For weak disorder, the frequencies of the vibrational modes follow the same 
dispersion law of a homogeneous chain, i.e., uj (k) = Wmax |sin(fi:/2)| with k 
being the mode wavenumber and Wmax = \/4:k/M the largest frequency. 

Eq. (!2T|) has the form of the stationary Schrodinger equation for the ID 
Anderson model. One can use this analogy to jump to the conclusion that, 
in an inflnite chain, the inverse localisation length of the vibrational modes 
is given by the expression 

i'\^) = 7^t^^'(^) Win) (22) 



2M2 
where 

oo 

Vr(K) = 1 + 2^x(0cos(2/k) (23) 
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is the power spectrum of the random succession {Snin} [3]. Making use of 
the dispersion law, the inverse locahsation length f l22p and the power spec- 
trum fl23|) can be expressed as functions of the frequency u rather than of the 
wavenumber k. As can be seen from eq. (I22l) . the inverse localisation length 
vanishes whenever the power spectrum W{k) is zero. Specific long-range 
correlations of the disorder correspond to the power spectra which vanish in 
continuous frequency intervals: this observation is the key to obtain delocal- 
isation transitions (within the second-order approximation) in ID disordered 
models [5|. 

The theoretical formula fl22|) for the inverse localisation length is valid 
only for an infinite chain, but the delocalisation transition which it predicts 
persists also in finite chains with N > 300 masses, as shown by numerical 
computations. This can be appreciated in Fig. [H which represents the in- 
verse entropic localisation length of the vibrational modes as a function of fre- 
quency. One can clearly see two windows of localised modes separated by the 
intervals of extended modes. The entropic localisation length, first applied 

Figure 1: Inverse entropic localisation length versus Wfe/wmax- The inset 
shows the entropic localisation length versus Uk/uj^aax- The data were 
obtained for a chain with = 500 masses and disorder strength a'^/M'^ = 
0.01. 
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in [T7] to quantum chaos, is a measure of the effective number of the eigen- 
vector components significantly different from zero. Specifically, the inverse 



entropic localisation length for an eigenmode with displacement vector e^"' is 



1 2 



log 



(i) 
2 



IS 



defined as lN{^i) = exp [SN{uJi)] where SN{(^i) - 
the Shannon entropy associated to the z-th mode. For localised eigenmodes, 
the entropic localisation length is proportional to the localisation length com- 
puted in the infinite-chain limit. Fig. [T]shows the inverse entropic localisation 
length for a random chain for which the binary correlator ( 120|) has the form 



x(0 



[sin (2K2O ~ ^i'^ (2ki/)] 



2 - Ki) I 

exhibiting the power-law decay typical of long-range correlated disorder. 



(24) 



4 Anomalous transport properties of the chain 
with correlated disorder 

We can now analyse how the delocalised modes affect the heat transport 
along the chain. The effect is best understood in terms of the modal heat 
flows ( fTSl) and (fT9l) . In figs. [2] and [3] we show the average modal flows in 
the stationary regime for fixed and free boundary conditions (the average 
was taken over 1000 disorder realisations). In both figures we show the 
modal fluxes for several kinds of harmonic chains: a homogeneous chain, a 
chain with uncorrelated disorder, and three chains with long-range correlated 
disorder. The latter cases of correlated disorder, which we refer to as case 1, 
2, and 3, have localised modes in two windows [ui : UJ2] and [u^ : 004], and 
extended modes in the complementary frequency intervals. In case 1, the 
limits of the localisation windows are cbi = 0.35, 0)2 = 0.4, 0)3 ~ 0.916, 0)4 ~ 
0.937, where a), = Wj/wmax- In case 2 we enlarge the windows of localised 
modes, setting their limits at cDi = 0.3, 0)2 = 0.5, 0)3 ~ 0.866, 0)4 ~ 0.954. 
Finally, in case 3 we merge the localised- mode intervals by setting a;2 = ^3; 
we thus generate a unique window [co'i,ci;4] of localised modes with edges 
ui = 0.3, UJ4 ^ 0.954. 

In all cases the chain consists of = 500 masses and the disorder strength 
is cT^/M^ = 0.01. As can be seen from figs. |2] and |3l the heat transport 
is effectively suppressed in the windows of localised modes, being close to 
that of the homogeneous chain in the windows of extended modes. By an 
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Figure 2: Average modal fluxes Jfc/Jmax versus normalised frequency Uk 
(fixed boundary conditions). Jmax is the largest modal flux. 
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Figure 3: Average modal fluxes J^/ Jmax versus normalised frequency (free 
boundary conditions). Jmax is the largest modal flux. 
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appropriate selection of the windows of localised and extended modes, one 
can either enhance (as in case 1 and 2) or suppress (as in case 3) the heat 
flow as compared with the case of uncorrelated disorder. 

The delocalisation of normal modes induced by long-range correlations 
of the isotopic disorder has repercussions on the conductivity of the random 
chain. The most striking consequence is that, by increasing the number of 
delocalised modes in the medium-to-high frequency regions, one can obtain 
random chains which behave more and more as homogenous chains. In par- 
ticular, the conductivity of such random chains scales as x oc N, regardless 
of the chosen boundary conditions. 

To understand this effect, let us focus on the case in which long-range 
correlated disorder produces extended modes in a frequency window [u2, OJ^]. 
In the limit iV ^ 1, the total heat flow can be approximated with the sum 
of the modal fluxes (fT8|) over the frequency window [c(;2, ws]. Assuming that 
the delocalised eigenmodes are not substantially altered with respect to the 
homogeneous case (a conjecture confirmed by numerical simulations for the 
weak-disorder case), one obtains that the heat flow in the stationary regime 
is 



\kB{T+ - T_) ^ r (,)] 2 \kB{T+ - T_) 



2M ^ 

k=k2 



2M 



(25) 



By dividing the flow f l25|) for the temperature gradient (T+ — T_)/N, one 
obtains that the effective conductivity of such a chain scales as ~ A^, as 
is the case for a homogeneous chain, and in contrast with the behaviour of a 
totally random chain. For the latter, heat transport is essentially due to the 
0(A^^/^) low-frequency modes and one has ~ N~^^'^ for fixed boundary 
conditions and ~ N^^"^ for free boundary conditions [T|. 

Our predictions are confirmed by the numerical data. To analyse the 
scaling law of the conductivity, we assumed an asymptotic behaviour of the 
form X = AN"' with a and A constants. We numerically computed the 
conductivity x for chains of length ranging from N = 200 to = 1000 
(with disorder strength cr'^/M'^ = 0.01). By fitting the data for logx versus 
log A^, we obtained values close to one for the exponent a. Specifically, in the 
case of fixed boundary conditions we respectively obtained = 0.93±0.01 



(2) 

and = 0.84 ± 0.02 for the cases 1 and 2 of correlated disorder. In 

„(i) 



the same cases, for free boundary conditions, the fitting values were a^^^^ 

(2) 
^free 



(2) 

0.96±0.01 and aL = 0.87±0.01. These results show that, when long-range 
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correlations delocalise a finite fraction of medium-to-high frequency modes 
(as in case 1 and 2), the enhanced conductivity does scale as in the case of a 
homogeneous chain. Especially in case 1, with a larger number of delocalised 
modes, the critical exponent takes values very close to unity and is essentially 
independent from the boundary conditions. 

We would like to stress that the enhancement or reduction of the conduc- 
tivity manifests itself in the time domain, through the change of the speed 
with which the system relaxes to equilibrium or approaches the steady-state 
regime. By increasing the number of delocalised modes in the medium-to- 
high frequency region, one can significantly speed up the approach to the 
steady-state condition. This can be seen in fig. H] which represents the evo- 
lution of the temperature of the central mass of a chain of = 501 atoms 
connected to two baths at the temperature T = 100. All the masses of the 
chains were initially at rest in their equilibrium positions, and the evolu- 
tion of the temperature was computed using eq. (|T6|) . The data shown in 




Fig. Hlwere obtained for the coupling constant A = 0.01 and disorder strength 
a^/M^ = 0.01. Fixed boundary conditions were chosen; the results, however, 
turn out to be qualitatively the same for free boundary conditions. 
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5 Conclusions 



This work analyses the heat conduction in a random chain coupled to ther- 
mal baths. We focus on the problem of how specific long-range correlations 
of the disorder influence the thermal properties of the model. We show that 
long-range correlations of the disorder can delocalise a finite fraction of vi- 
brational modes in the middle region of the frequency spectrum even for 
chains of relatively modest length (A^ > 300). These delocaliscd modes mod- 
ify both the approach to the stationary regime and the steady-state thermal 
properties of the chain. The larger the number of extended modes, the faster 
the chain tends to reach the steady-state condition and the closer the ther- 
mal properties are to those of a homogeneous chain. In particular, when the 
conductivity is enhanced it exhibits anomalous scaling x oc iV" with a ~ 1, 
independently of the boundary conditions. The long-range correlations of the 
disorder can also be used to strengthen the localisation of the chain eigen- 
modes: in this way one obtains a chain which is a better thermal insulator 
than a random chain with uncorrelated isotopic disorder. 

l.F.H.-G. and L.T. gratefully acknowledge the support of the CONACyT 
grant No. 84604 and of the CIC-2009 grant (Universidad Michoacana) . The 
work of F.M.I, was partly supported by CONACyT grant No. 80715. 
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